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Spallation neutron production in proton induced reactions on Al, Fe, Zr, W, Pb and Th targets 
at 1.2 GeV and on Fe and Pb at 0.8, and 1.6 GeV measured at the SATURNE accelerator in Saclay 
, is reported. The experimental double-differential cross-sections are compared with calculations 

performed with different intra-nuclear cascade models implemented in high energy transport codes. 
The broad angular coverage also allowed the determination of average neutron multiplicities above 
, 2 MeV. Deficiencies in some of the models commonly used for applications are pointed out. 
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PACS numbers: PACS numbers: 25.40.Sc, 24.10.Lx, 25.60.Dz, 29.25.Dz 



I. INTRODUCTION 



Large numbers of neutrons can be produced through spallation reactions induced by an intermediate energy (around 
1 GeV) proton accelerator on a heavy element target. With the progress in high intensity accelerators it is now possible 
to conceive spallation sources that could compete with high flux reactors. Several spallation sources for solid state and 
material physics are under construction or study in the USA (SNS Q), in Europe (SINQ 0], ESS |3|) and in Japan 
^ 1 (NSP Q). Spallation neutrons can also be used in Accelerator Driven Systems (ADS) to drive sub-critical reactors, 
in which long-lived nuclear waste could be burnt [|5], [| or energy produced (tJ. All these systems have in common 
a spallation target made of a heavy, either solid (W, Pb) or liquid (Hg, Pb, Pb-Bi eutectic) metal in a container 
r"| , (generally steel) which is separated from the vacuum of the accelerator by a thin window. 

A detailed engineering design of a spallation target needs a precise optimisation of its performances in terms of 
useful neutron production and a proper assessment of specific problems likely to occur in such systems, like induced 
radioactivity, radiation damage in target, window or structure materials, additional required shielding due to the 
presence of high energy neutrons, etc... This can be done by using Monte-Carlo transport codes describing the 
interaction and transport of all the particles created in nuclear reactions occurring inside the system to be designed. 
Generally, a high energy transport code (often based on the HETC code from ||), in which elementary interactions 
are generated by nuclear physics models, is coupled below 20 MeV to a neutron transport code like MCNP || that 
uses evaluated data files. Although the spallation mechanism has been known for many years, the models used in such 
codes, Intra-Nuclear Cascade followed by evaporation-fission, have never been really validated on experimental data 
and large discrepancies remain both between experimental data and model predictions and between different models. 



This was particurlarly obvious from the OECD/NEA intercomparisons (lfj, 11, [L2| of these codes, regarding neutron 
and residue production. This led to the conclusion that many improvements of the models are still needed but also 
that there was a lack of experimental data to make a good validation, especially above 800 MeV. Among the needed 
data, the energy and angular distributions of spallation produced neutrons are essential for model probing: the high 
energy part of the spectrum allows the testing of the intra-nuclear cascade while the low energy part of the spectrum 
is sensitive to the excitation energy at the end of the intra-nuclear cascade stage and the evaporation model. They 
are also important to optimize the target geometry since secondary particles contribute to the propagation of the 
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inter-nuclear cascade in a thick target and high energy neutrons are responsible for radiation damage in target and 
structural materials. 

During the last years, a wide effort has been made in several laboratories to measure spallation data reg arding 
neutron multiplicity distributions |lj| |l4], [T^J [lfj, light charged particles jlTj and heavy residues JT^, |o[ ||l[] in 
order to establish a base for the test and validation of the spallation physics models. We report in this paper on 
neutron production double-differential cross-sections, measured at the SATURNE synchrotron, induced by 0.8, 1.2 
and 1.6 GeV protons impinging on different targets. The experimental setup, although already discussed in detail 
in |2^, |23|1 , is presented in section 1 while the results are displayed in section 2. The aim was to measure at one 
energy (1.2 GeV) neutron spectra on nuclei representative of different parts of the periodic table of elements and at 
the same time corresponding to materials used in targets or structures of accelerator driven systems: Al, Fe, Zr, W, 
Pb and Th. At the two other energies only Fe and Pb targets were studied. Angular distributions covering 0° to 
160° were obtained and allowed the determination of average neutron multiplicities per reaction above 2 MeV. The 
Pb data have already been published in a letter j24| but are again reported here for the sake of completeness. Section 

3 is devoted to the comparison of the data to different intra-nuclear cascade models. 

II. EXPERIMENTAL APPARATUS 

The slow extraction of the beam delivered by the Saturne synchrotron did not allow conventional time-of-flight 
measurement using the HF signal of the accelerator. Therefore, the time-of-flight had to be measured between the 
incident proton passing through a thin scintillator placed in the beam and the detected neutron. For the highest 
neutron velocities this method becomes highly imprecise due to the limited available flight path, and a measurement 
of the (n,p) charge exchange by a spectrometer technique was adopted. These two methods had previously been 
tested in a first experiment in which neutrons were measured at 0° and have been respectively described in detail 
in two previous papers |2^, [2^] . Then, a new experimental area was designed to allow complete angular distribution 
measurements. The scheme of this experimental area is shown in Fig. [lj The beam comes from the left, hits the studied 
target and is deflected to a beam stopper (composed of lead and tungsten blocks) by a dipole magnet. To ensure the 
detection of the neutrons emitted by the target and not those coming from the beam stop, a large shielding of heavy 
concrete was built around the target. Because the most energetic neutrons are emitted in the forward direction this 
shielding is thicker between 0° and 90° (3.5 meters) than at backward angles (2.5 meters). The shielding is pierced 
by 12 circular holes at 0° and every 15° from 10° to 160°. Each channel is composed of two consecutive cylinders 
of respectively 110 mm and 136 mm in diameter with a length of respectively 1600 mm and 1800 mm at forward 
angles and 1100 mm and 1400 mm at backward angles. The solid angle is determined by the size of the neutron 
detectors, not the collimator. The two detection methods were used in dedicated runs due to the fact that they 
required different beam intensities: between 10 9 and 10 11 particles per second for the spectrometer (because of the 
low neutron-to-proton charge exchange efficiency) while the time-of-flight was limited to less than 10 6 particles per 
second by the use of the in-beam scintillator. The spectrometer was not used at angles larger than 85° since very 
high energy neutrons are expected only at forward angles. 

A. The time-of-flight method. 

The Saturne synchrotron was delivering continuous beam during spills of ~ 500ms with a repetition frequency 
depending on the beam energy (typically 1.5 seconds for 800 MeV protons in this experiment). A classic time-of-flight 
measurement between the neutron detector and a HF signal of the accelerator was thus not possible. We just recall 
here the principle of the method and the modifications performed since the first experiment at 0° [ p2[ . 

The time-of-flight is measured between the incident proton, tagged by a thin plastic scintillator (SC), and a neutron 
sensitive NE213 liquid scintillator (see Fig ||). The beam intensity is fixed at a maximum of 10 6 particles per second 
so that individual incident protons could be counted by SC. The target-detector distance depends on the angle but is 
about 8.5 meters. Up to ten angles can be explored simultaneously using neutron detectors composed of a cylindrical 
cell of NE213 liquid scintillator coupled to a photomultiplier. Six of them are cells of the multi-detector DEMON |2[| 
and the other four (called DENSE) are smaller detectors, optimized for low energy measurements. The latter detectors 
allow energy measurements with a reasonable precision from 2 to 14 MeV, while the DEMON cells are used between 

4 and 400 MeV. The characteristics of the DEMON and DENSE detectors are given in table [j| 



The DEMON detectors are placed in a shielding of paraffin loaded with Borax and Lithium to reduce the background. 
The DENSE detectors, smaller and thus less sensitive, do not need such a protection. The energy threshold of the 
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FIG. 1: Experimental area with time-of-flight and spectrometer setup. 




STOP 



START 



FIG. 2: Principle of the time-of-flight method shown here for the 0° measurement with the dipole magnet ensuring the deflection 
of the incident beam. 



detectors is adjusted using the Compton edge of the gamma rays delivered by a 22 Na and a 137 Cs radioactive source. 
The detection thresholds of the DENSE and DEMON are tuned to 1.0 and 1.9 MeV respectively. This allows 
measurements with a sufficiently well defined efficiency above 2 and 4 MeV respectively. 

These detectors are sensitive to neutrons, 7-rays and charged particles. A plastic scintillator NE102 placed in front 
of each counter (see Fig ||) tags events induced by a charged particle. The neutron-gamma discrimination is performed 
by a pulse shape analysis done as follows p5| : the charge delivered by the photomultiplier is measured by a QDC 
1612F during two different time intervals, a prompt one (125 ns long) and a delayed one (185 ns long delayed by 
65 ns) giving two charge values (Of) and (Qs) respectively. A bidimensional spectra Qf vs Qs allows the separation 
of neutrons and gammas (see Fig. [|). 

The neutron detector starts a gate of 500 ns, longer than the time-of-flight of 2 MeV neutrons (440 ns to cover the 
8.5 meters from the target to the detector). The signal from the scintillator (SC) in the beam is delayed by 500 ns 
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Characteristics 


DENSE 


DEMON 


Liquid scintillator 


NE213 


NE213 


Cell diameter 


127 mm 


160 mm 


Cell length 


51 mm 


200 mm 


Photomultiplier type 


9390 KB 


XP 4512 


Detector threshold 


1. MeV 


1.9 MeV 



TABLE I: Neutron detector characteristics. 




FIG. 3: Neutron gamma discrimination. 



and is the stop of the time measurement. Due to the beam intensity, although limited to less than 10 6 p/s, several 
protons can be detected during the gate, but only one has induced the detected neutron. For a common start (the 
neutron) up to ten stop signals are converted and stored with a multistop module (LECROY TDC 3377). The sum 
of these time spectra contains the real and accidental events. The accidental contribution is determined by the time 
measurement of unmatched start and stop signals. The background effect is taken into account by a measurement 
with an empty frame at the location of the target. 

The knowledge of the neutron detector efficiency being crucial for this experiment, measurements and calculations 
have been performed to determine it over the whole energy range. From 2 to 17 MeV, measurements were made at the 
Bruyeres-le-Chatel Van de Graff accelerators as described in p2| . The quasi mono-energetic neutrons are produced 
by 7 Li(p,n) 7 Be, 3 H(p,n) 3 Hc, 3 H(d,n) 4 He and 2 H(d,n) 3 He reactions and the efficiency is determined by comparison 
with a standard detector (full triangles in fig.[| and fig. ||). 

At higher energies (30 < E < 100 MeV), experiments have been performed at the TSL Uppsala facilities in 
Sweden |2(|. A neutron beam is produced by 7 Li(p,n) 7 Be reaction in the [100-180 MeV] range. The neutron detector 
efficiency is then measured using n-p elastic scattering and the simultaneous detection of the correlated n and p; the 
so called associated particle method. For energies from 150 to 800 MeV, the d + Be break-up reaction is used at 
Saturne to produce quasi monoenergetic neutrons. The deuteron beam intensity is measured by activation of a carbon 
foil [^7j and the neutron flux is deduced from d + Be break-up cross-sections Q . 

The results are displayed in Fig. |i| and Fig. |[ For the Saturne measurement, the neutron flux attenuation along 
the 7 to 8.5 meters of flight path in air and in the veto scintillator is the same during the efficiency measurements 
and the real experiment: therefore it has not to be corrected for. On the contrary, experimental points measured 
at Bruyeres and Uppsala have to be corrected to take into account the difference in neutron flux attenuation due to 
the difference in distance and the absence of the veto detector. The corrected efficiencies used in the data analysis 
are represented by the solid line in Fig. |] and Fig. || for DEMON and DENSE detectors respectively. The efficiency 
calculations performed with 05S |^9| (diamonds) and a modified version of the KSU (triangles) Monte-Carlo codes, 
agree fully with experimental results. In the original KSU code |]30| some reactions are not taken into account and 
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NE213 threshold 1.9 MeV 
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Neutron energy (MeV) 

FIG. 4: Efficiency of a DEMON neutron detector as a function of energy measured at Bruyeres-le-Chatel, Uppsala and Saturnc 
(black symbols) and calculated with 05S [^9| and a modified version of the KSU §(J code (grey symbols). The solid line 
represent the final parameterization used in the data analysis after correction of the neutron flux attenuation in air along the 
flight path in Saturne and in the NE102 scintillateur in front of each detector. 



40.0 




Neutron energy(MeV) 



FIG. 5: Efficiency of a DENSE neutron detector as a function of energy. The symbols represent the values measured at 
Bruyeres-le-Chatel and the dashed line a fit of these values. The solid line is the final parameterization used in the data 
analysis (after correction of the neutron flux attenuations in Saturne set-up). 

at 90 MeV the sum of the cross sections over all the inelastic processes gives 165 mb, that is 90 mb less than the 
global cross section measured by Kellog |pl] . In our modified KSU code, the total inelastic cross-section has been 
normalised to the Kellog measurements by an appropriate increase of the (n, a) light response. This is justified by 
the fact that the missing reactions in the KSU model produce essentially recoil nuclei like deuterons, tritons, alpha, 
lithium or boron whose light production is close to that of alphas. The line presented in Fig. ^ is the efficiency used for 
all DEMON detectors. The DENSE detectors being used only between 2 and 14 MeV, efficiency has been determined 
by measurements only in the 2-17 MeV range (see Fig. ||). The efficiency was assumed to be the same for all the 
detectors of the same type. 

The beam is monitored by the start scintillator located in front of the target. Uncertainties on the cross-section 
determination are due to statistical and systematic errors. Systematic errors come mainly from the subtraction of 
accidental coincidences, gamma rejection and efficiency determination. The error on the latter mainly depends on 
the absolute calibration procedure used at the different accelerators and are estimated to be 10% for the Bruyeres 
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measurement (knowledge of the standard detector), 4% at Uppsala (n — p cross-section) and 10% at Saturne (d + Be 
cross-section). The values are summarized in table 0. 



Error origin 


2-20 MeV 


20-100 MeV 


100-400 MeV 


subtraction of accidental coincidences 


5.8% 


5.8% 


5.8% 


gamma rejection 


2.9% 


2.9% 


2.9% 


efficiency determination 


10% 


4% 


10% 


TOTAL uncertainty 


11.9% 


7.6% 


11.9% 



TABLE II: Estimations of systematic errors in the time-of-flight method. 



The neutron energy resolution depends on a time and a geometrical component and is given by : 

^ = 7(7+1) 

with 7 the Lorentz factor 
the time resolution 
the geometrical component 

The time resolution (a t — 1.5 ns and 2.0 ns for DEMON and DENSE detectors respectively) is estimated by the 
measurement of the FWHM of the prompt gamma peak on the time-of-flight spectra. The geometrical component 
comes from the target thickness (1 to 3 cm) and from the size of the sensitive area of the detector. The interaction 
probability being constant as a function of the depth the standard uncertainty is er; = Q . Thus 07 = 6 cm and 
1.5 cm for DEMON and DENSE detectors respectively. 
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FIG. 6: Geometrical and time contributions to the energy resolution as a function of neutron energy for a DEMON detector. 
The time-of-flight length is 8.5 m, the time uncertainty 1.5 ns and the length uncertainty is 6 cm. 

The energy resolution is plotted as a function of the energy in figure It appears clearly that above 400 MeV, 
this time-of-flight method doesn't allow neutron energy measurement with a resolution better than 12%. Therefore 
another complementary detection system has been developed and is described in the next section. 

B. Proton recoil spectrometer. 

High energy (i.e. above 200 MeV) neutrons are detected using (n,p) scattering on a liquid hydrogen converter and 
detecting the emitted proton in a magnetic spectrometer (Fig mj . We present here only a comprehensive description 
of the measurement with a special emphasis on the modifications compared to the the first experiment at 0° which is 
detailed in ref Q . 

The spectrometer is composed of the dipole magnet VENUS, which generates a 0.4 T field, and of 3 multiwire 
proportional chambers, Ci, C2 and C3, of respective active area 20x20cm 2 , 80x40cm 2 and 100x80cm 2 . Each chamber 
is composed of 2 sets of wires equipped with PCOS2 electronics allowing the localization in the X — Y plane. The 
wires are spaced by 1.27mm for Ci and 2mm for C2 and C3. 
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The acquisition is triggered by the coincidence between the plastic scintillator SI and the large wall of NE102 plastic 
scintillators behind VENUS. This wall is made of 20 horizontal slats with a photomultiplier on each side. 




FIG. 7: Principle of the high energy neutron measurement (E > 200 MeV) with the spectrometer, shown here at 0°. 

A second dipole magnet, CHALUT, deviates in the horizontal plane the charged particles created in the target, in 
the concrete or in air. The field integral of this magnet was OATm during our experiment. The thin plastic scintillator 
SA tags the possible charged particles remaining in front of the liquid hydrogen target. The spectrometer (hydrogen 
target and detectors) and CHALUT are placed on platforms which could rotate from to 85° only since, for larger 
angles, very few neutrons with energies higher than 400 MeV are expected. A single setting of the magnetic field in 
VENUS (0.4T) was chosen during the standard measurements. The center of the chamber C3 was shifted by 40 cm 
to the left with respect to the beam axis in order to optimize the detection of deflected protons. 

The trajectory reconstruction of the charged particles emitted by the hydrogen target is deduced from the impact 
coordinates in Ci, C2 and C3. The well known magnetic field inside VENUS gives their momentum. The geometric 
calibration of the multiwire chambers was done using a low intensity 800 MeV proton beam without the hydrogen 
target and successive magnetic fields of —0.2, 0., 0.2 and 0.4 Tesla. 

The liquid hydrogen target is a cylinder of 12.8 cm in diameter. It has a useful thickness of 0.94 g/cm 2 and is 
located 8.45 meters from the production target, covering an angular aperture of 0.43°. The entrance and outgoing 
windows are titanium foils, 100 [im thick . 

Several types of charged particles are created in the H2 target through the following processes : np— >np, np— »np7r°, 
np— >pp7r - , np— >d7r°, np— >2n7r + . The charged particles with different masses are identified using the biparametric 
representation of the ST Wall time-of-flight versus the momentum measured with the wire chambers and VENUS 
(Fig. |). 

The incident proton beam is monitored by two telescopes viewing a 50 fj,m mylar foil placed upstream in the beam 
(~ 20m from the target). The absolute calibration of these telescopes is obtained by a comparison with the activation 
of a carbon sample |27j] measured in a dedicated run. A calibration with activation of Al foils was also done and gives 
a very comparable result. 

The response function of the spectrometer which takes into account the contribution of elastic and inelastic processes 
arising in the hydrogen target has been measured with quasi-monoenergetic neutrons produced by the break-up of 
deuterons or 3 He beams on a Be target. From various beam energies, neutrons between 0.2 and 1.6 GeV were 
produced for this calibration. The neutron flux is obtained from the known (n,p) elastic scattering cross-sections |33[ ] 
and the normalized response functions are then used to unfold the measured proton spectra. This procedure gives the 
normalized neutron energy distribution. It is described in detail in ref [ p3| . The maximum energy available at Saturne 
(1.6 GeV/ A for 3 He) and the growing importance of the inelastic processes set a limit to this unfolding procedure. 



Systematic errors in this method arise mainly from the beam calibration, the spectrometer response function and 
the unfolding procedure. The estimations of these 3 errors are given in table III. They are less than 10% at 800 
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Momentum (GeV/c) 



FIG. 8: Bidimensional spectra to identify protons, deuterons and pions produced in the hydrogen target. 



Error origin 


Beam Energy 


0.8 and 1.2 GeV 


1.6 GeV 


Beam monitoring 


< 5.8% 


< 5.8% 


Spectrometer response function 


< 4% 


< 11.5% 


Unfolding procedure 


< 5.8% 


< 8.6% 


TOTAL uncertainty 


< 9.1% 


< 15.5% 



TABLE III: Estimations of systematic errors as a function of incident energy. 



and 1200 MeV but reach 15.5% at 1600 MeV due to the increasing inelastic contribution. Error bars on the results 
presented in this paper take into account only statistical uncertainties except at the very low energies where the 
increase of uncertainty associated with the proximity of the detection threshold is added. 



C. Targets 

The same targets were used for both the time-of-flight and spectrometer methods. Because of the low beam intensity 
imposed by the detection of the incident proton with the SC scintillator, we had to use rather thick targets in order 
to keep a significant counting rate. They were 3cm diameter cylinders made of natural material with thickness in the 



centimeter range, shown in Table IV, depending on the elements 



Target 


Al 


Fe 


Zr 


W 


Pb 


Th 


Thickness (cm) 


3 


3 


3 


1 


2 


2 



TABLE IV: Target thicknesses for the different materials used in this experiment. 
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FIG. 9: Results obtained at 1600 MeV on the iron target with 
the DEMON detectors and the spectrometer in the overlap- 
ping region at 10 and 25°. 



FIG. 10: Id. but for the overlapping region of the DENSE 
and DEMON detectors at close angles 40-55° and 145-160° 
since evaporation neutrons are emitted nearly isotropically. 



III. DOUBLE-DIFFERENTIAL CROSS-SECTIONS 



A. Experimental results 



In order to show the consistency of the three different sets of detectors, figs. ^ and [l(J display details of the 
double-differential cross-sections obtained for a Fe target at 1600 MeV. In fig. || data obtained at 10 and 25° with the 
DEMON detectors (filled circles) and the spectrometer (squares) are shown. It can be seen that in the overlap regions, 
i.e. between 200 and 400 MeV, the data are compatible within the error bars. Actually, 1600 MeV corresponds to the 
worst case since, as mentioned above, the spectrometer unfolding procedure is approaching its limits of reliability. 
Other examples of the good agreement between both methods regarding Pb at 800 and 1200 MeV were shown in [f24[ . 
In fact, for all the measurements, the data from the DEMON detector and spectrometer always agree within less 
than 15% at 1600 MeV and 10% at lower energies. Concerning the comparison of spectra obtained between 4 and 12 
MeV with the DENSE and DEMON detectors, measurements at the same angles were done only for Pb at 1200 MeV 
and also shown in ]24| . However, since below 15 MeV neutrons mostly come from an evaporation process which is 
practically isotropic in the laboratory system, results obtained at near angles can be compared. This has been done 
in fig. |l(] where 25 and 145° DEMON spectra are plotted together with respectively 40 and 160° DENSE ones. Both 
sets of data are consistent and agree within less than 10%. This appears to be verified whatever the target and the 
energy. Therefore, a single set of data merging the different measurements by taking their mean values in the overlap 
regions has been processed and will be shown in the following. 

Fig. [ll] shows comparisons of our data with previously obtained ones by Amian et al. p4j and Nakamoto et al. |3f| 
using time-of- flight. In each case we compare the data at the closest possible angles, adding when appropriate 
results from DEMON and DENSE detectors at two different angles. Actually, our 10, 25-40, 55, 85-100, 115-130 and 
145-160° data are displayed together with the previous ones at respectively 15, 30, 60, 90, 120 and 150°. Our Pb 
measurements at 800 MeV (left) fully agree with Amian ones, as already noticed in |24|], while Nakamoto cross-sections 
are systematically lower at low neutron energy and higher at energies between 10 and 100 MeV. For Fe (center), only 
data from (3^] are available and we observe a slightly less good agreement between the two works. However, it should 
be stressed that contrary to |Q our (and ref |35|) targets are not really thin and secondary reactions increase the 
number of low energy neutrons (as discussed below). This effect is visible only below 4 MeV and appears to be larger 
for a 3cm thick Fe target than for a 2cm thick Pb one, as shown by the simulations in fig. [l2| This could explain 
why we measure more neutrons than Amian et al. at low energies in the case of iron. At 1600 MeV for Pb (fig. |ll| 
right), we have compared our data to ref pj| data obtained at 1500 MeV. This is possible since, from our 1200 and 
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FIG. 11: Comparison of our results with the data from Amian et al. j34j and Nakamoto et al. |55| at different close angles; 
left : for Pb at 800 MeV; center : for Fe at 800 MeV; right : for Pb at 1600 MeV. Each successive curve, starting from the 
smallest angle, is scaled by a multiplicative factor of 10 _1 . 

1600 MeV measurements, we could infer that cross-sections should differ by less than 5% between 1500 and 1600 
MeV, apart from the high energy part of the spectra at very forward angles. As observed at 800 MeV, we get higher 
cross-sections below 7 MeV and lower ones at intermediate energies. Since the thickness of the targets are the same in 
both cases, this can be understood only by differences in the neutron detector efficiency determination. As mentioned 
earlier, at intermediate energies, our experimentally determined efficiency [^6| is higher than the one calculated using 
the standard KSU code and thus we used a modified version. As far as we know, it was the standard version of KSU 
that was used to determine the detector efficiency in j3fj| ]. This could explain the discrepancy. At high energies, with 
the spectrometer, we obtain a much better energy resolution than in ref. [[35] (because of their limited flight path) 
that allows us to distinguish structures due to direct reactions at forward angles. 

As already mentioned, the thickness of our targets induces some distortion in the double-differential spectra because 
of the slowing down of the incident proton and the probability that some of the energetic emitted particles may undergo 
secondary reactions. The first point was discussed in |24| and results in a slight shift and a broadening of the quasi- 
elastic and inelastic peaks at very forward angles. The second effect is expected to lead to a depopulation of the 
high and intermediate energy parts of the spectra and an increase of the number of low energy emitted neutrons. 
Calculations using the LAHET high-energy transport code system |3(| (using Bertini as intra-nuclear cascade model 
and pre-equilibrium) were performed for both a target with the actual geometry and an infinitely thin one in order to 
assess the order of magnitude of the effect. In fig. [jj], results are shown for Pb and Fe at 800 MeV. It can be seen that 
the difference is very small for the 2cm lead target, for which it is significant only between 2 and 3 MeV, and a little 
larger for the 3cm iron one. In the latter case, the disappearance of intermediate energy neutrons is also perceptible. 
Similar results are found with the other targets and at other energies, the effect being maximum for Fe and W. 

All the measured angular distributions are presented together with the model calculations in the next section. Data 
have been taken at 0.8, 1.2 and 1.6 GeV on Pb and Fe targets and at 1.2 GeV on Th, Pb, W, Zr, Fe and Al targets 
at 0, 10, 25, 55, 85, 130 and 160°. 

B. Comparison with models 

Spallation reactions are generally described by a two step mechanism: a first stage in which successive hard 
collisions between the incident particle and the individual nucleons of the target nucleus lead to the emission of a few 
fast nucleons, then, the decay of the excited remnant nucleus by emission of low energy particles or, sometimes for 
heavy nuclei, by fission. The first step is generally described by Intra-Nuclear Cascade models while evaporation-fission 
models are used for the second one. Some authors introduce a pre-equilibrium stage between intra-nuclear cascade and 
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Energy (MeV) Energy (MeV) 



FIG. 12: Effect of the target thickness at 800 MeV for Pb (left) and Fe (right) : calculations done with the 
LAHET code ]36j for an infinitely thin target (solid line) and a 2cm Pb and 3cm Fe target respectively (dashed 
line) . 



de-excitation. In high energy transport codes, the most widely used intra-nuclear cascade model is the old Bertini 37] 
one dating from 1963. However, several other models are available, such as the Isabel ]38| and the Cugnon 4C] 
INCL models, which have brought some improvements in the physics. The most widely used eva por ation model in 
the domain of spallation reactions is the Dresner model |4lJ , usually associated with the Atchison |4^| fission model. 

The high energy part of the neutron spectra enables one to directly probe the intra-nuclear cascade models. Low 
energy neutrons, which are the majority of the neutrons produced in spallation reactions, are emitted during the 
evaporation process. However, their number mainly depends upon the intra-nuclear cascade stage since the cascade 
determines the initial excitation energy of the decaying hot residue and, therefore, the number of evaporated particles. 
Actually, evaporation neutron spectra for a given excitation energy are not expected to depend very much on the 
evaporation-fission model contrary to light charged particle spectra or to residual nuclei production for which emission 
barriers and competition between the different decay modes are not so well established. This is why we have made 
comparisons with different intra-nuclear cascade models in order to test their validity and understand their differences, 
using always the same evaporation-fission model. 

All the calculations discussed in the following have been done with high energy transport codes in which the actual 
thickness and diameter have been taken into account. In order to have sufficient statistics, calculations were done 
for angular bins of 5 degrees (except at degrees where it is only 2.5 degrees). The INCL model does not predict 
a correct total reaction cross-section mainly because the diffuscness of the nuclear surface is not taken into account. 
Therefore, the INCL calculations were renormalized to the total reaction cross-sections given by the Bertini model 
which appear to be in very good agreement with experimental values from |48|. 

In |24j], for lead, we presented calculations performed with the TIERCE |43| high-energy transport code system 
developed at Bruyeres-le-Chatel (which is very similar to LAHET) using either the Bertini or the Cugnon INCL model 
with the same evaporation- fission model (based on the Dresner- Atchison model). It was shown that, at the three 
measured energies, the Bertini model was largely overpredicting the experimental data while INCL was giving a rather 
good agreement. This was ascribed to the higher excitation energy, E*, obtained at the end of the cascade stage with 
the Bertini calculation than with INCL. This assumption can be verified in fig. [l3| where the E* distribution obtained 
with Bertini (solid line) is shown to extend to much higher values than INCL (dashed-dotted line) and gives also a 
higher average value (265 versus 213 MeV). These calculations were performed for thin targets at 1 GeV. The same 
observations were also made in [|~7| and |4£] where a similar plot was shown for p+Au reactions and INCL was found 
to give the best agreement with the excitation energy deduced from the neutron multiplicity distributions. Several 
reasons can explain the difference in E* between the two models : first, INCL leads to the emission of more pions 
than Bertini. However, the difference in average E* due to the energy carried away by the pions is only 30 MeV. 
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FIG. 13: Excitation energy distribution in the p (1 GeV) + 
Pb reaction found with the Bertini (solid line), Bertini+pre- 
equilibrium (dashed line), Isabel (dotted line) or INCL 
(dashed- dotted line) intra- nuclear cascade models. 



FIG. 14: Excitation energy distribution in the p (1 GeV) + 
Pb reaction calculated with the Cugnon INCL model using 
the standard (solid line) or a strict Pauli blocking (dashed 
line). 



Second, as mentioned in p4[ , the Pauli blocking is treated in a different way. In Bertini, only collisions of nuclcons 
with momentum larger than the Fermi momentum are allowed while, in INCL, the actual phase space occupation 
rate is taken into account. This leads to a less stringent condition, therefore more cascade particles can escape and 
make the energy remaining in the nucleus lower. This is illustrated in fig. |lj where the E* distributions with a strict 
Pauli blocking (as in Bertini) and the standard one are shown. The decrease is obvious. 

In the LAHET code system |36| it is possible to add after the intra-nuclear cascade stage a pre-equilibrium Q| 
stage which is expected to reduce the excitation energy of the nucleus by emission of intermediate energy particles 
prior to the evaporation. Besides, this is also the recommended option by the LAHET authors |Q. Also available is 
the Isabel model which can be used only up to 1 GeV. As can be seen in fig. [l3[ both models lead to excitation energy 
distributions close to the one found with INCL. Isabel is used with the partial Pauli blocking (recommended) option, 
which, as in INCL, is supposed to take into account the depletion of the phase space due to the emission of cascade 
particles. Here, we show calculations performed with both models and the same Dresner-Atchison evaporation-fission 
at 800 MeV for the Pb and Fe target. In the following, Bertini plus pre-equilibrium will be referred to as BPQ. 
Fig. [l^ presents the calculated neutron spectra compared to the experimental data. It can be observed that, for Pb, 
the BPQ calculation reproduces very well the data, except at very forward angles and high neutron energies where 
the peak corresponding to the excitation of the A resonance appears much too high. This is a deficiency of the 
Bertini intra-nuclear cascade model, already pointed out in 46 as due to a bad parameterisation of the NN — > NA 
reaction angular distribution. The problem does not exist with Isabel. Both models correctly predict the low energy 
part of the spectra. This can be understood by the respective excitation energy distributions being similar in their 
extension to the one found with INCL (see fig. [l3]). The high energy neutrons above 85° are also well reproduced by 
both calculations but Isabel underestimates cross-sections at backward angles in the intermediate energy region. For 
iron, Isabel presents the same features as for lead while BPQ now overpredicts low and intermediate energy neutron 
production at forward angles, indicating that the angular distribution of pre-cquilibrium neutrons is probably too 
much forward-peaked. 

At 1200 MeV, the use of Isabel in LAHET being limited to 1 GeV, we compare the data with only BPQ and 
INCL calculations in figs. |l6|, |l7| and [l8[ for all the targets. We also performed calculations, which are not shown 
here, using Bertini without pre-quilibrium. Whatever the target, this model yields too many low energy neutrons, 
emphasizing that it leads to too high excitation energies. For Th, Pb and W, both BPQ and INCL models give 
a reasonable agreement with the data, although BPQ tends to slightly overestimate the production of intermediate 
energies neutrons. As the target becomes lighter, this trend is amplified and BPQ begins to also overpredict low 
energy cross-sections. This is an indication that the addition of a pre-equilibrium stage after intra-nuclear cascade to 
decrease the too large excitation energy found in Bertini may not be the proper solution: in fact, it seems difficult to 
obtain the correct evaporation neutron production without overestimating intermediate energy cross-sections (which 
are enhanced by pre-equilibrium). On the contrary, INCL reproduces quite well the results for all the targets, proving 
that the model has a correct mass dependence. Only for the light targets at very backward angles, the high energy 
neutron production is underpredicted. 

At 1600 MeV, fig. |l| displays the results for the Pb and Fe targets. For BPQ the tendencies noticed at 1200 MeV 
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FIG. 15: Experimental p (800 MeV) + Pb (left) and Fe (right) neutron double- differential cross-sections 
compared with calculations performed with LAHET using either Bertini plus pre-equilibrium (solid line) or 
Isabel (dashed line) intra-nuclear cascade model. Each successive curve, starting from 0°, is scaled by a 
multiplicative factor of 10 _1 . 

are growing worse: even for Pb, the agreement is not very good between 10 and 40 MeV. Since the high energy part 
of the spectra is always rather well reproduced (except at 0°), this seems to point out a wrong dependence of the 
pre-equilibrium emission also with incident energy. Here again, INCL gives a satisfactory agreement with the data 
for both targets. 

In summary, we can conclude that INCL is able to globally reproduce the bulk of our data, with some slight 
discrepancies in the angular distributions. The Bertini model followed by pre-equilibrium, although it is found to be 
an improvement compared to Bertini alone, works well for Pb at 800 MeV but fails as the energy is increased and the 
target gets lighter. 



IV. AVERAGE MULTIPLICITIES PER REACTION NEUTRON 



Since our double-differential cross-sections nearly cover the full angular range with sufficiently close measurements, 
it has been possible to infer from the data average neutron multiplicities per reaction above the energy threshold of our 
detectors. This has been done by interpolating between the measured angles, integrati ng ov er 4-7T and then dividing by 
the reaction cross-section taken from [Eg) and the result is shown in tables |\], [Vj and VII for 2 different energy bins, 
corresponding roughly to evaporation and cascade neutrons respectively. Since we divide by the reaction cross-section, 
the multiplicities obtained are numbers of neutrons per primary reaction and therefore contain the effect of secondary 
reactions. The interpolation between angles was done assuming that the angular dependence of the cross-sections 
was the same as the one calculated by the TIERCE code. This was necessary in particular between 2 and 4 MeV 
where we have only a few points from the DENSE detectors. The uncertainty on these interpolations was assessed 
by using different intra-nuclear cascade models in TIERCE and different interpolation procedures. The values given 
in the tables take into account this uncertainty plus the systematic errors discussed in section 2. The experimental 
values are compared with the neutron multiplicities given by the two codes, TIERCE-INCL and LAHET-BPQ. For 
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FIG. 16: Experimental p (1200 MeV) + Th (left) and Pb (right) neutron double- differential cross-sections 
compared with calculations performed with LAHET using either Bertini plus pre-equilibrium (solid line) or 
INCL (dashed line) intra-nuclear cascade model. 

the calculations, 0-2 MeV and total multiplicities are also given. Results are also shown in fig. ^ as a function of 
the mass of the targets studied at 1200 MeV and in fig. ^l] versus incident energy for the iron and lead targets. The 
comparison of the experimental to calculated average neutron multiplicities confirms in a rather concise way what 
has been observed in the preceeding section. In all cases, INCL agrees with the data within the error bar while BPQ 
tends to overpredict 2-20 MeV neutron multiplicities, i.e. evaporation neutron production, especially at 1200 and 
1600 MeV. For high energy neutrons (above 20 MeV) the sensitivity to the models is less important. This is most 
likely because of compensating effects, BPQ predicting more intermediate energy neutrons because of pre-equilibrium 
while INCL spectra often extend to higher energies. However, a significant deviation from the experiment is found at 
the two highest energies for iron with BPQ. 

Neutron multiplicity distributions have been measured by the NESSI collaboration ^?], using the technique of 
a gadolinium loaded liquid scintillator tank. The efficiency of this type of detector is large for low energy neutrons 
and decreases rapidly above 20 MeV. Mean neutron multiplicities were obtained by NESSI at 1200 MeV on W and 
Pb targets with the same thickness as ours. Although the threshold of our detectors did not allow us to measure 
neutrons with energies lower than 2 MeV, it is tempting to compare the two results. For Pb at 1200 MeV, the NESSI 
collaboration measured for a 2cm thick, 15cm diameter target, 14.6 neutrons per reaction which after efficiency 
correction j|7| |^)| amount to 20.3. In it is indicated that the average multiplicities have not been corrected for 
the additional neutrons coming from secundary reactions in the liquid scintillator. This effect was however investigated 
and estimated to be 5% for the 2cm thick Pb target (see fig. 8 of |49]). Taking this correction into account, the NESSI 
multiplicity is then 19.3. In our case, we can estimate the total number of neutrons by adding the experimental values 
measured between 2 and 1200 MeV to the [0-2 MeV] multiplicity given by the codes. If we take the average between 
INCL and BPQ values we find 16.9 neutrons per reaction. The error is estimated to be of the order of 15% taking 
into account the errors discussed above plus the extrapolation to low energies. While the thickness of the targets in 
both experiments is the same the diameter is somewhat larger in the NESSI case. We performed a simulation with 
LAHET to investigate this effect. It appeared that because of secondary reactions of particles emitted sidewards the 
number of neutrons is 5% larger than with a 3cm diameter target. This means that for a 15cm diameter target we 
would find 17.8 ± 2.7 neutrons to be compared to 19.3 ± 1.9 if we take an uncertainty of 10% for the NESSI result. 
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FIG. 17: Experimental p (1200 MeV) + W (left) and Zr (right) neutron double-differential cross-sections 
compared with calculations performed with LAHET using either Bertini plus pre-equilibrium (solid line) or 
INCL (dashed line) intra-nuclear cascade model. 



The same can be done for the 1cm thick W target. In this case, the additional contribution due to the target diameter 
is 3%. Therefore, we obtain 16.9 ± 2.5 to be compared to 18.0 ± 1.0 for NESSI [|^J after efficiency correction and if 
we assume, as in Pb, 5% additional neutron from secundary reactions in the liquid scintillator. As it can be seen in 
fig. ^2| we can conclude that our r esults are compatible with those from NESSI within the error bars. 

Also shown in tables [v|, [vj], VII are the averaged kinetic energies carried away by the neutrons, extracted from the 
double-differential cross-sections multiplied by energy using the same procedure as for the multiplicities and compared 
with the calculations. For the 2-20 MeV bin, conclusions similar to what was stated for multiplicities can be drawn, 
reflecting the fact that ExM„ is governed by M n in so far as the same evaporation model is used in both calculations 
and thus gives an identical energy spectrum for the low energy neutrons. For the high energy bin, the compensating 
effect noticed for the multiplicities seems to be even stronger and, regarding our uncertainties, it is not possible to 
discriminate between the two models. It is interesting, nevertheless, to remark that these high energy neutrons carry 
out the major part (from 80% for Th to 98% for Al) of the emitted neutron energy and a large amount (about 30%) of 
the incident proton energy In a thick target this will play an important role in the spatial distribution of the energy 
deposition and particle production. 



V. CONCLUSION 



In this paper, we have displayed double-differential cross-sections measured on a wide set of targets and at different 
energies. This has allowed a comprehensive comparison with some of the high energy models commonly used in 
high energy transport codes for applications. In particular, we have compared different models describing the first 
stage of the reaction (intranuclear cascade possibly followed by pre-equilibrium) keeping the same model for the de- 
excitation process. As we had shown in a previous paper that the Bertini intra-nuclear cascade was leading to too 
high excitation energies at the end of the first reaction stage, we have tried to use the option, available in LAHET 
and recommended by its authors, of adding a pre-equilibium stage. This largely improves the predictions of the code. 
However, discrepancies tend to appear and grow larger as the energy is increased above 800 MeV and as the target 



16 




FIG. 18: Experimental p (1200 MeV) + Fe (left) and Al (right) neutron double-differential cross-sections 
compared with calculations performed with LAHET using either Bertini plus pre-equilibrium (solid line) or 
INCL (dashed line) intra-nuclear cascade model. 
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TABLE V: Average neutron multiplicities per primary reaction and kinetic energy carried out by the neutrons (MeV) obtained 
by integration of the double-differential cross-sections at 800 MeV and compared with calculations using TIERCE-Cugnon or 
LAHET-BPQ codes. 



becomes lighter. The Isabel model was also tried at 800 MeV. It gave reasonable agreement with the lead data but 
less good one for iron. Unfortunately, since the use of Isabel is limited to 1 GeV in LAHET, it was not possible to 
test the energy dependence of this model. Finally, we have shown that the use of the Cugnon intranuclear cascade 
model, INCL, implemented in the TIERCE code, is able to fairly reproduce the whole bulk of our results. However, 
it should be recalled that this model still suffers from serious deficiencies mostly due to the fact this it does not treat 
correctly the diffuseness of the nuclear surface. This was the reason why we had to renormalize the calculation to 
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FIG. 19: Id. but for p (1600 MeV) + Pb (left) and Fe (right) 
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FIG. 20: Average neutron multiplicities per primary reaction at 1200 MeV for the different targets. Left: [2-20 
MeV] neutron multiplicities; Right: [20 MeV-Einc] neutron multiplicities. 



the correct total reaction cross-section. Also, the sharp surface approximation makes it impossible to have a correct 
prediction of the most peripheral collisions: this was clearly seen, for instance, in the isotopic distributions of residual 
nuclei close to the projectile presented in ref. fl9fl . A new version of the Cugnon model is in progress f3(if which is 
expected to solve this problem. 

All the data presented will be given to the EXFOR data base or are available on request. 
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FIG. 21: Average neutron multiplicities per primary reaction for Pb and Fe as a function of incident energy. Left: 
[2-20 MeV] neutron multiplicities; Right: [20 MeV-Einc] neutron multiplicities. 
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FIG. 22: Total neutron multiplicities per primary reaction at 1200 MeV for 15cm diameter W (1cm thick) and Pb (2cm thick 
targets obtained by the NESSI collaboration after efficiency and secundary reactions in the liquid scintillator corrections [^, ^ 
compared to the values estimated from our measurement. 
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TABLE VI: Average neutron multiplicities per primary reaction and kinetic energy carried out by the neutrons (MeV) obtained 
by integration of the double-differential cross-sections at 1200 MeV and compared with calculations using TIERCE-Cugnon or 
LAHET-BPQ codes. 



and Nucl. Sci. Eng. 49, 110 (1972). 
[9] J.F.Briesmeister, ed., MCNP A General Monte Carlo N-particle Transport Code, LA-12625-M, Los Alamos Nat. Lab. 
(1997). 

[10] M.Blann et al., International Code Comparison for Intermediate Energy Nuclear Data, OECD/NEA (1994). 
[11] D.Filges et al., OECD thick target benchmark, Report NSC/DOC(95)2 (1995). 

[12] R.Michel and P. Nagel, Int. Codes and Model Intercomparison for Intermediate Energy Activation Yields, OECD/NEA 
(1997). 

[13] L. Pienkowski et al, Phys. Rev. C 56, 1909 (1997). 
[14] X. Lcdoux et al, Phys. Rev. C 57, 2375 (1998). 

[15] B. Lott et al, Nucl. Instr. and Meth. in Phys. Res. A 414, 117 (1998). 
[16] D. Hilscher et al, Nucl. Instr. and Meth. in Phys. Res. A 414, 100 (1998). 
[17] M. Enke et al., Nucl. Phys. A 657, 319 (1999). 

[18] R. Michel et al, Nucl. Instr. and Meth. in Phys. Res. B 103, 183 (1995). 
[19] W. Wlazlo et al, Phys. Rev. Lett. 84, 5736 (2000). 
[20] T. Enqvist et al., Nucl. Phys. A 686 (2001) 481. 

[21] F. Rejmund et al., Nucl. Phys. A 683 (2001) 540; J. Benlliure et al., Nucl. Phys. A 683 (2001) 513. 



20 



Energy M° xp M* NCL Mg PQ || E x M* xp E x Mjf CL E x Mg PQ 

E inc = 1600 MeV 



Pb(p,xn)X <t r = 1717 mb 



0-2 MeV 


6.0 


6.6 




6. 


7. 


2-20 MeV 10.1 ±1.0 


9.9 


12.2 


65. ± 7. 


61. 


81. 


2Q-E inc . 3.4 ± 0.5 


3.1 


3.4 


410. ± 41. 


422. 


416. 


Total 


19.0 


22.2 




489. 


503. 


Fe(p,xn)X 








<JR = 


774 mb 


0-2 MeV 


1.2 


1.6 




1. 


1. 


2-20 MeV 1.9 ± 0.2 


2.3 


3.1 


14. ± 1. 


16. 


24. 


2Q-E inc . 1.9 ± 0.3 


2.0 


2.3 


341. ± 34. 


363. 


387. 


Total 


5.5 


7.0 




380. 


412. 



TABLE VII: Average neutron multiplicities per primary reaction and kinetic energy carried out by the neutrons (MeV) obtained 
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